The Glass Transition and Liquid-Gas Spinodal Boundaries of Metastable Liquids 
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A liquid can exist under conditions of thermodynamic stability or metastability within bound- 
aries defined by the liquid-gas spinodal and the glass transition line. The relationship between these 
boundaries has been investigated previously using computer simulations, the energy landscape for- 
malism, and simplified model calculations. We calculate these stability boundaries semi-analytically 
for a model glass forming liquid, employing accurate liquid state theory and a first-principles ap- 
proach to the glass transition. These boundaries intersect at a finite temperature, consistent with 
previous simulation-based studies. 
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When a substance is brought below its equilibrium 
freezing temperature and yet maintained in the liquid 
state, it is referred to as a supercooled liquid. Even if 
crystallization is avoided a liquid cannot be supercooled 
indefinitely; it transforms into an amorphous solid via 
the glass transition, at a temperature determined in ex- 
periments by the cooling rate. The ideal glass transition 
refers to the putative transformation to an amorphous 
solid which has been argued to occur in the limit of in- 
finitesimally small cooling rates. Such a transition would 
define an ultimate limit beyond which a substance cannot 
exist in the supercooled liquid state. 

If a liquid were instead to be heated beyond the boil- 
ing point - and the transformation into the gas phase 
averted - it is termed a superheated liquid. Such a liquid 
cannot be superheated to arbitrarily high temperatures 
since growing density fluctuations will trigger a catas- 
trophic transformation into the gaseous state at rates 
which increase with the degree of superheating. Mean- 
field theory identifies a degree of superheating for which 
the compressibility diverges, at which point the liquid 
would spontaneously transform into a gas. The locus of 
such points defines another ultimate limit to the stability 
of the liquid state, the liquid-gas spinodal p). 

Clarifying the relationship between these two limits 
is fundamental to a deeper understanding of the liquid 
state. In addition to conventional materials in the liq- 
uid state, much interest has recently focussed on such 
questions in the context of colloidal fluids. In these sys- 
tems, a unified picture of arrested states (colloidal gels 
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and glasses) has been the subject of considerable investi- 
gation |2| . The understanding of the role of the unstable 
region in the phase diagram demarcated by the liquid- 
gas and gas-liquid spinodals is a necessary component of 
a comprehensive picture of these arrested states Q ■ 

In a previous investigation Q of the model liquid stud- 
ied in this letter, through computer simulations and the 
energy landscape approach, it was found that the glass 
transition line and the liquid-gas spinodal intersect at a 
finite temperature, thus predicting a glass-gas limit of 
mechanical stability at lower temperatures. More re- 
cently, similar results have been obtained from simula- 
tions of a model of ortho-terphenyl0 as well as calcu- 
lations based on model energy landscapes 0, Q- In this 
letter, we describe an analysis of these limiting lines for a 
realistic glass former, which uses (i) accurate liquid state 
theoretical methods to calculate the equation of state 
of the system, and thus the liquid-gas spinodal, coupled 
with (ii) a first-principles approach to the glass transi- 
tion proposed by Mezard and Parisi 0, 0] to evaluate 
the glass transition line. Our calculations are quantita- 
tively in reasonable agreement with the previous study, 
and also predict that the glass transition line and the 
liquid-gas spinodal intersect at a finite temperature. Our 
results thus provide strong support for the proposal that 
such an intersection of these limiting lines is the generic 
behavior for liquids. 

The model system we study is the Kob-Andersen bi- 
nary mixture Lennard- Jones (KA BMLJ) liquidj^l, con- 
sisting of an 80:20 mixture of A and B particles in- 
teracting via a Lennard-Jonnes potential with parame- 
ters e AB /e AA = 1.5, e BB /e AA = 0.5, cr AB /a AA = 0.8, 
o~bb/o- aa = 0.88. We report energies in units of e AA , 
lengths in units of o~ AA , and temperatures in units of 
e AA /k B . This system has been extensively studied as a 
model glass former. Our procedure for calculating the 
glass transition line, which we summarize below, is close 
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to that employed in Ref.s [l3l.ll4| . 

We compute the equation of state of the system we 
study, KA BMLJ, from evaluation of the A - A, A - B 
and B — B pair correlation functions. Liquid state the- 
oretical methods provide accurate analytic tools for the 
computation of such pair correlation functions and we ex- 
ploit such methods here. Our calculations use the Zerah- 
Hansen (ZH) closure scheme |9j coupled to the Ornstein- 
Zernike (OZ) equation. The Zerah-Hansen closure inter- 
polates between hyper-netted chain (HNC) and the soft- 
core mean-spherical approximation (SMSA) closures [loj: 
as described in Ref. 13], such an interpolation reproduces 
thermodynamic quantities (energy and pressure) accu- 
rately. 

For the binary system, with indices a and b indicating 
particle types, the OZ equation is 



h a b(k) = Cab(k) 



E 



PiCai(k) * hib(k) 



(1) 



where pi are the partial number densities, h a b is de- 
fined as g a b{r) — 1 (g a b(x) being the ab pair correlation 
function) and c a b{f) is the direct correlation function. 
h a b{k) etc indicate space Fourier transforms. Defining 

lab{r) = h ab (r) - c a b(r), 
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ab(k) = p t Cai(k)(cib(k) +^ib{k)). 



(2) 
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The ZH closure is of the form 



c a c»p 2 g ab (r) 



exp(-/3<(r)) 



(3) 



exp[f ab (r)[ lab (r) - Pv a A b (r)} -1 
f ab (r) 



where (3 = 1/fcsT, c a and c b are the concentrations of 
A type and B type particles, respectively, vr and va are 
attractive and repulsive contributions of the pairwise po- 
tential v, and are given by vr = v(r) — v mm (r < r min ) 
and v B , = (r > r min ); v A = v mm (r < r min ) and 
Vr = v(r) (r > r m i n ). The function / provides an 
interpolation between the SMSA and the HNC: f ab = 



1 



exp 



with a governing the switch between 

HNC and SMSA. The parameter a is typically chosen 
by demanding thermodynamic consistency between equa- 
tions of states obtained by the virial and compressibility 
routes. Other choices, however, are possible an( i 
here we choose the parameter a by comparison with the 
simulation data for the studied system in a range of tem- 
perature and density values. We find, as expected, that 
the value of a that yields the best match of energy and 
pressure values with those from simulations varies both 
with temperature T and density p. Based on the estima- 
tion of such a values for a small set of temperature and 
density values, we choose a linear functional form for the 
(3 and p dependence of a, a(p, T) = 11.17 — 11.7/5+6.65/3, 
which allows us to obtain quantitatively good results for 
thermodynamic quantities. 



We solve the OZ equation, along with the ZH closure, 
for t a h using the procedure outlined in Ref. [l5j. Pair 
correlation functions are obtained first for low density 
and high temperatures. The pair correlation functions so 
obtained are used as initial guesses for the iterative solu- 
tion of the OZ and ZH equations at lower temperatures 
and higher densities. We require iterative convergence at 
the level of 1(T 9 in J2i(9n+i (n) - gn(n)) 2 , with n being 
the iteration number and i the real space mesh index. 

In Fig. 1, we compare isotherms obtained using 
this scheme for temperatures T — 0.4, 0.5, 0.6 with the 
isotherms derived, in Q, from an empirical equation of 
state based on simulation data. Data shown in Fig. 1 
demonstrate that the equation of state from the present 
calculations are in good quantitative agreement with the 
results in Figure 2 shows the potential energies (e pot ) 
obtained from the present calculations as a function of T 
and p, along with available values from computer simu- 
lations [4|, demonstrating excellent agreement with the 
simulation results. The comparisons in Fig.s 1 and 2 in- 
dicate that the pressure and energy values obtained from 
the present calculations are quantitatively quite accurate. 

In Fig. 3, we show the isotherms obtained for temper- 
atures in the range T = 0.4 — 0.9. We obtain the location 
of the spinodal at each temperature from the condition 

9P _ n 
dp - u - 
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FIG. 1: Pressure vs density isotherms obtained using the 
Zerah-Hansen scheme (symbols), compared with the equation 
of state obtained from fits to the simulation data (lines) for 
temperatures T — 0.4, 0.5 and 0.6. 

The glass transition temperatures at each density are 
calculated u sing approach proposed by Mezard and Parisi 
[Til IT2I liH Il7| . which has been employed to calculate 
the glass transition line for the KA BMLJ in [l3| . In this 
calculation, assuming the existence of an exponentially 
large number of free energy minima whose degeneracy 
determines the configurational entropy or complexity, the 
thermodynamics of a system of m copies of the system 
which are restricted to be in the same free energy min- 
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where Su q is the entropy of the liquid, and 



Sglass((3) 



3N 



l+log(2ic))--<Tr(log(/3M))>, (5) 



where M is the Hessian matrix of the form M/^^u^ 



(n -rjfe) - v iiv {r l - tj) with v p 



The 



indices /i and v run over the coordinates and i , j and fc 
run over particle indices. The trace is calculated using 
the harmonic resummation scheme, details of which (and 
of the discussion above) can be found in 

PEam. as 

in 0, , the ideal glass transition Tk is located with 
the condition S con f{TK) = 0. 



FIG. 2: Comparison of potential energy obtained using the 
Zerah-Hansen scheme (lines) with data from molecular dy- 
namics simulation, shown for densities p — 1.1, 1.15, 1.2, 1.25 
as a function of the inverse temperature /3. 
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FIG. 4: Configurational entropy S con f for densities p = 
1.0, 1.05, 1.10, 1.15, 1.20 are shown. For low temperatures the 
Sconf is obtained via an extrapolation shown by the dashed 
lines. The Kauzmann temperature is identified with the tem- 
perature at which the configurational entropy vanishes. 



FIG. 3: Pressure vs. density isotherms obtained using the 
Zerah-Hansen scheme for temperatures varying from T = 0.4 
to T = 0.9. The spinodal lines are obtained from the location 
of the minima along these isotherms. 



imum is considered. With the replicas remaining con- 
fined close to each other, the free energy can be written 
in terms of their center of mass coordinates and the ex- 
cursions of individual replicas from the center of mass. 
Evaluation of the contribution from the latter is done 
within approximation schemes pertaining to the Hessian 
or potential energy second derivative matrix. The con- 
figurational entropy is obtained by differentiating the re- 
sulting free energy per replica with respect to their num- 
ber, in the limit of m — > 1. The resulting expression for 
the configurational entropy is 



The entropy of the liquid is obtained from the internal 
energy (e(/3) = + e pot ) via 

f fi • ■ 

S liq ({3) = S? iq +0e Hq - / d(3e(f3). (6) 
Jo 

Here S° iq is the entropy of the perfect gas in the infinite 
temperature limit and is given by S° iq = 1 — log(p) — 
c a log(c a ) — c b log(c b ) where c a is the concentration of A 
particles and c b is the concentration of B type particles. 

Both Suq and S g i ass require knowledge of the liquid 
state pair correlation functions. However, at low tem- 
peratures the liquid integral equations do not converge, 
and hence extrapolations (of the form aT -2 / 5 + b for Su q 
and a + b log(T) for S g i aS s, based respectively on theo- 
retical predictions for dense liquids |l9j | and the harmonic 
approximation to the glass free energy) are used to obtain 
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FIG. 5: Limits of stability of the liquid phase from present 
work - spinodal labeled T s and glass transition labeled Tk 
- compared with previous work |4|: (circles) spinodal esti- 
mated from inverse compressibility (dashed line is a guide to 
the eye), (solid line) from the empirical equation of state, and 
(long dashed line) glass transition line from energy landscape 
approach. The present calculation of the spinodal agrees quite 
well with that from the empirical equation of state, while the 
glass transition temperatures are seen to be slightly overes- 
timated. The two limiting lines intersect at p ~ 0.95, as 
compared to p ~ 1.08 from earlier work. (The intersection 
density in Q is 1.02 if the spinodal estimate from the empir- 
ical equation of state is used instead of the simulation esti- 
mates, a number closer to the one obtained from our present 
calculations.) 



Sconf at low temperatures. Figure 4 shows the configu- 
rational entropy vs. T for p — 1.0,1.05,1.1,1.15,1.20. 



The liquid-gas spinodal and the glass transition lines 
obtained are shown in Fig. 5, along with results from ear- 
lier work Q . Our results reproduce the features observed 
in earlier work both qualitatively and, to a good extent, 
quantitatively, thus providing support for the general ap- 
plicability of this approach. The calculated spinodal tem- 
peratures are in reasonable agreement with the ones es- 
timated earlier from the empirical free energy. The glass 
transition temperatures obtained here are slightly higher 
than those estimated in previous work; this may be a 
consequence of the approximations involved in evaluat- 
ing the glass entropies in the present case. 

Together with previous work 0, H, El 1 these results 
lend strong support to the scenario that the liquid-gas 
spinodal and the glass transition line should typically in- 
tersect at a finite low temperature. Unlike previous work, 
where the results were based either on simulations, or 
on assumed model energy landscapes, our calculations 
describe a realistic model glass former, whose proper- 
ties are evaluated within a self-contained and accurate 
framework. We believe that the approach used here and 
its extensions may be useful in understanding arrested 
states in contexts similar to the one we study. 
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